#include <stdexcept>
#include "prot-filer/xtc_filer.h"

namespace prot_filer
{
using namespace mili;

enum Axis
{
    XAxis,
    YAxis,
    ZAxis,
    // -------------------
    AxisSize  // dejar este valor a lo ultimo
};

typedef float EjesCartesianos[AxisSize];

struct AtomWrapper
{
    AtomWrapper(BasicProtein& prot) :
        protein(prot),
        atom_number(0)
    {};
    inline static float scale(float coord)
    {
        return 10.0f * coord;
    }
    inline void set_atom(rvec atom)
    {
        protein[atom_number].x = scale(atom[0]);
        protein[atom_number].y = scale(atom[1]);
        protein[atom_number].z = scale(atom[2]);
        protein[atom_number].vdw = (atom_number % 3 == 0) ? N : ((atom_number % 3 == 1) ? CA : C);   // Check this!!!
        atom_number++;
    }
    BasicProtein& protein;
    unsigned int atom_number;
};

struct ProteinWrapper
{
    ProteinWrapper(Protein& prot) :
        protein(prot)
    {};
    inline static float scale(float coord)
    {
        return coord;
    }
    inline void set_atom(rvec atom)
    {
        const float x = scale(atom[0]);
        const float y = scale(atom[1]);
        const float z = scale(atom[2]);
        protein.addItem(Coord3d(x, y, z));
    }
    Protein& protein;
};

static const float _precision = 1000.0;
static const float box[DIM* DIM] =
{
    9.0, 0.0, 0.0,
    0.0, 9.0, 0.0,
    0.0, 0.0, 9.0
};

void XtcWriter::initialize()
{
    for (size_t i(0); i < DIM; ++i)
    {
        for (size_t j(0) ; j < DIM; ++j)
        {
            default_box[i][j] = box[i * DIM + j];
        }
    }
}

bool XtcWriter::_open(const std::string& name)
{
    xfp = xdrfile_open(name.c_str(), "w");
    step_number = 1;
    return xfp != NULL;
}

void XtcWriter::take_step()
{
    ++step_number;
    _time += _TIME_STEP;
}

void XtcWriter::writextc(BasicProtein& patm, int nres)
{
    const size_t natom = ATOMS_PER_RESIDUE * nres;
    auto_array<EjesCartesianos> cxtc(new EjesCartesianos[natom]); /* cosas para escribir el xtc*/

    for (size_t i = 0; i < natom; i++)
    {
        cxtc[i][XAxis] = patm[i].x / 10.0;
        cxtc[i][YAxis] = patm[i].y / 10.0;
        cxtc[i][ZAxis] = patm[i].z / 10.0;
    }
    write_xtc(xfp, natom, step_number, float(step_number), default_box, cxtc.get(), 10000.0);
}

void XtcWriter::write(BasicProtein& patm, const AnglesData& angles_data)
{
    checkOpen();
    writextc(patm, angles_data.nres);
    take_step();
}

void XtcWriter::write(const Protein& protein)
{
    checkOpen();
    // Taken from xtc_filer.cpp, part of project parallel-clusterer http://parallel-clusterer.googlecode.com.
    write_xtc(xfp, int(protein.items()), step_number, _time , default_box,
              const_cast<rvec*>(reinterpret_cast<const rvec*>(&protein.front())) , _precision);
    // const_cast is safe here, after analysis.
    // The base address of the first Coord3d in the proteins works as the rvec vector.
    take_step();
}

void XtcWriter::write(const std::vector<Coord3d>& backbone)
{
    write_xtc(xfp, int(backbone.size()), step_number, _time , default_box,
              const_cast<rvec*>(reinterpret_cast<const rvec*>(&backbone.front())), _precision);
    // const_cast is safe here, after analysis.
    //Idem, &backbone[0] works fine
    take_step();
}

void XtcWriter::_close()
{
    if (xfp != NULL)
    {
        xdrfile_close(xfp);
        xfp = NULL;
    }
}

bool XtcReader::_open(const std::string& name)
{
    xfp = xdrfile_open(name.c_str(), "r");
    filename = name;
    auto_array<char> fname(new char[filename.size() + 1]);
    strcpy(fname.get(), filename.c_str());
    int result = read_xtc_natoms(fname.get(), &_atoms_in_a_protein);

    if (exdrOK != result)
        throw std::runtime_error("Error initializing database. Wrong filename?");

    if (NULL == xfp)
        throw std::runtime_error("Error opening file.");

    return xfp != NULL;
}

bool XtcReader::read(BasicProtein& protein)
{
    checkOpen();
    AtomWrapper atom_wrapper(protein);
    return read(atom_wrapper);
}

bool XtcReader::read(Protein& protein)
{
    checkOpen();
    ProteinWrapper protein_wrapper(protein);
    return read(protein_wrapper);
}

template<class T>
bool XtcReader::read(T& protein)
{
    auto_array<char> fname(new char[filename.size() + 1]);
    strcpy(fname.get(), filename.c_str());
    read_xtc_natoms(fname.get(), &_atoms_in_a_protein);
    auto_array<rvec> _atoms_vector(new rvec[_atoms_in_a_protein]);
    int    step = 0; // To avoid valgrind error.
    float  time;
    float  prec;
    Matrix box;
    // read one protein
    int result = read_xtc(xfp, _atoms_in_a_protein, &step, &time, box, _atoms_vector.get(), &prec);
    bool end = false;
    if (exdrENDOFFILE != result)
    {
        if (exdrOK != result)
        {
            throw std::runtime_error("ERROR: read_xtc " + result);
        }
        else
        {
            for (int atom_number(0); atom_number < _atoms_in_a_protein; ++atom_number)
            {
                protein.set_atom(_atoms_vector[atom_number]);
            }
        }
    }
    else
    {
        end = true;
    }
    return end;
}

void XtcReader::_close()
{
    if (xfp != NULL)
    {
        xdrfile_close(xfp);
        xfp = NULL;
    }
}
}
